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Abstract. - We investigate stall force and polymerization kinetics of rigid protofilaments in a 
microtubule or interacting filaments in bundles under an external load force in the framework 
of a discrete growth model. We introduce the concecpt of polymerization cycles to describe the 
stochastic growth kinetics, which allows us to derive an exact expression for the stall force. We 
find that the stall force is independent of ensemble geometry and load distribution. Furthermore, 
the stall force is proportional to the number of filaments and increases linearly with the strength 
of lateral filament interactions. These results are corroborated by simulations, which also show a 
strong influence of ensemble geometry on growth kinetics below the stall force. 



Introduction. — Polymerization of cytoskeletal fila- 
ments is essential for various cellular processes, such as 
motility or the formation of cellular protrusions includ- 
ing filopodia or lamellipodia 1,2 . Single polymerizing 
filaments can generate forces in the piconewton range, as 
has been demonstrated experimentally for microtubules 
(MTs) [3] . Such force generation mainly relies on the gain 
in chemical bonding energy upon monomer attachment 
[3] . An opposing force slows down filament growth, which 
finally stops at the stall force representing the maximal 
polymerization force a filament can generate. Therefore, 
the stall force is the essential quantity to characterize poly- 
merization forces. 

Cellular force generating structures such as filopodia are 
made of polymerizing ensembles of interacting actin fila- 
ments [S]. Particularly important are bundles of parallel 
filaments, which can hold together by crosslinking pro- 
teins or unspecific attractive interactions. Stall forces of 
polymerizing actin bundles could be determined experi- 
mentally only recently 6 . 

MTs are tubular filaments, which also consist of an en- 
semble of typically 13 interacting protofilaments (PFs). 
The force velocity-relation of polymerizing MTs has been 
experimentally determined in refs. [Hll^, where stall forces 
around 5 pN have been obtained. 

An ensemble of many non-interacting filaments or PFs 
is believed to have higher stall forces than a single fil- 
ament because of load sharing effects. First fits of the 
experimental data on MT growth in ref. [3] were based on 
the assumption of load sharing and application of ratchet 



models for a single rigid PF 8. . An explicit continuous 
model for N rigid PFs in a MT under load resulted in 
stall forces cx N^^"^ 9 . For an analogous discrete growth 
model of PFs in a MT it has been shown that the stall 
force of N PFs increases oc N compared to a single PF 
[To] , in agreement with equal load sharing. Variants of 
this model which allow a better fit of experimental data 
were discussed in |llj . 

In addition to load sharing effects, crosslinking or at- 
tractive lateral interactions within filament bundles or be- 
tween PFs can allow the ensemble to generate even higher 
forces by exploiting the additional interaction energy |12] . 
For fiexible filaments, zipping mechanisms for force gener- 
ation can even rely exclusively on attractive interactions 
[l3] . These results suggest that the stall force of interact- 
ing filaments or PFs increases by the additional interaction 
energy per length that a bundle gains upon assembly. Lat- 
eral interactions between PFs of a polymerizing MT have 
been considered in refs. [T4U16) . Also for MT growth, ap- 
proximative analytical results in refs. [14l[15] suggest that 
the stall force of interacting PFs increases by the addi- 
tional interaction energy per length that the MT gains 
upon assembly. 

Apart from this progress, an exact result for the stall 
force could not be derived so far. Furthermore, the geom- 
etry of the bundle or tube, i.e. the relative positioning of 
filaments or PFs in the ensemble, has an impact on the 
mechanics of monomer insertion under load and on the 
lateral interactions, which are involved. 

In this Letter, we investigate the combined effects of at- 
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tractive filaments interactions and ensemble architecture 
on the growth kinetics under a compressive force using a 
discrete growth model [TOllMllTS] . Based on the concept of 
polymerization cycles we derive an exact analytical result 
for the stall force. This result shows that for the discrete 
model introduced in refs. [TnillHIin] the stall force is a 
universal quantity, which only depends on the polymer- 
ization energy gain and the interaction strength between 
filaments or PFs. The stall force is independent of ensem- 
ble geometry and independent of the distribution of load 
force and interaction energy between attachment and de- 
tachment rates. The result also shows that the stall force 
increases oc N, i.e., linearly in the number of filaments. 

Using stochastic simulations based on the Gillespie al- 
gorithm we find that the growth kinetics below the stall 
force depends sensitively on the strength of interactions 
between rigid filaments in a bundle or PFs in a MT and 
on the geometry of the bundle or tube. We find differ- 
ent shapes of the force-velocity relation as well as a com- 
plex non-monotonic dependence of the growth velocity on 
the relative filament positioning. This dependence is very 
pronounced at low forces well below the stall force but 
vanishes upon approaching the stall force resulting in a 
geometry independent stall force. 

Our results are of particular interest with respect to the 
growth kinetics of MTs, which usually contain 13 PFs. 
Our results imply that a two-start helical structure (with a 
helical pitch of one tubulin dimer) , which is often assumed 
in modelling, has a distinct force-velocity relation but an 
identical stall force as the actual three-start helix structure 
(with a helical pitch of three tubulin monomers) found by 
electron microscopy [17) . 

Model and simulation. — We consider a filament 
consisting of N rigid PFs in a tube-like arrangement such 
that each PF has two neighbors and periodic boundary 
conditions apply. For actin bundles, this model neglects 
effects from thermal shape fluctuations [T^ and the exis- 
tence of defects within the bundle structure p"9H2Tj . Each 
PF consists of monomers of size d, see fig. [1] The attach- 
ment and detachment rates for monomers are fcon and fcoff , 
respectively, and related to the polymerization energy gain 
Ep > upon adding a monomer by fcon/fcoff = gEp/keT 
at temperature T. For MTs, each monomer is a tubulin 
dimer, and we will neglect hydrolysis of GTP such that fcon 
and fcoff are attachment and detachment rates for GTP- 
tubulin dimers. We also neglect catastrophes and consider 
MTs only in their growing phase. Effects of hydrolysis are 
shortly discussed in the end. 

Each PF has attractive lateral interactions with its 
two neighbors; the corresponding lateral association en- 
ergy (per length) is e; > 0. Thus, apart from the poly- 
merization energy, an attaching monomer gains the lat- 
eral interaction energy £;A£, where A£ > is the addi- 
tional contact length with neighboring monomers, which 
is created upon inserting the monomer. Likewise, a de- 
taching monomer looses a corresponding interaction en- 



ergy resulting in A£ < 0. If the PF can equilibrate 
its remaining configurational degrees of freedom suffi- 
ciently fast during addition or removal of a monomer, 
thermodynamics requires that lateral interactions change 
on- and off-rates of that monomer such that fc*,-j/fc*g = 
(fco„/fcoff)exp(£,|A£|/fcBr). 

The influence of an external load is described by a force 
F, which acts only on the leading PF. We denote the po- 
sition of the tip of the leading filament by x. Insertion of 
a monomer changes the position of the tip of the leading 
PF by Ax > in the on-process, which gives rise to an 
additional mechanical energy F\x. Likewise, removal of 
a monomer in the off-process gives rise to Ax < 0. If we 
assume again that the PF can equilibrate its remaining 
configurational degrees of freedom sufficiently fast during 
addition or removal of a monomer under force, thermody- 
namics requires further modification of on- and off-rates 
such that 

Kn_ _ {Ei\Al\-F\i:^x\)/kBT t-,^ 

toff '^oS 

If on- and off-rates fco„ and k*^ are specified separately 
the thermodynamic constraint ([ij allows to introduce a 
load distribution factor 9p and a lateral energy distribu- 
tion factor di, 

h* - k p(-eFF\Ax\ + (l-ei)ei\Al\)/kBT 

fc*g = fc^^g((l-e.)F|Ax|-9,e,|A£|)/fc«T_ (3) 

In general, load and energy distribution factors can de- 
pend on the specifics of insertion and removal and differ 
for each polymerization step. A reasonable assumption is 
that the external load only affects the on-rate (correspond- 
ing io Op — 1) and that the on-rate is diffusion-limited 
and not affected by lateral interactions (corresponding to 
9i = 1) but other choices are thermodynamically possible. 
In the absence of lateral interactions, this model (with 
9p — 1) was introduced by van Doorn et al. [10]. Lateral 
interactions have been included in refs. [HHTB] . 

The geometry of the bundle or tube has an impact on 
the mechanics of monomer insertion under load and on 
the lateral interactions, which are involved. Therefore, we 
expect a strong influence of geometry on growth kinetics. 
We control the ensemble geometry by shifting the relative 
position of neighboring PFs in the initial configuration by 
a distance h, which we call a geometry parameter, see fig. 
d] Because of the periodic boundary conditions, we have 

— 1 relative displacements of h between neighboring PFs 
and one relative displacement of (N — l)h. For /i = we 
have an aligned or "flat" initial configuration with all PF 
tips at the same height. For symmetry reasons, geometry 
parameters h and —h are equivalent. The parameter h can 
also be shifted by multiple monomer sizes d corresponding 
to the insertion of additional monomers without changing 
the kinetics. Therefore, it is sufficient to consider values 
< h < d/2. In the following we will measure h in units 
of d and use h = h/d. 
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Fig. 1: Schematic picture of the model of interacting PFs 
polymerizing against the external load force F for A'' = 4. Po- 
sitions of monomers of size d in neighboring PFs are shifted 
by the offset h. Insertion of a monomer (yellow) might change 
the position of the leading PF by Ax (right picture) and create 
additional contact length with neighboring PFs (red). 



Of special interest are MTs, which are usually built from 
iV = 13 PFs with tubulin dimers of size c? ~ 8 nm J/7j. 
Often it is assumed that MTs with = 13 PFs have 
an offset h = 1/13 [IHIIIII, which results in a symmetric 
arrangement without a seam in the structure, where the 
tube closes. The actual structure exhibits a three-start 
helix with a helical pitch of three tubulin monomers (or 
1.5d) resulting in h = 1.5/13 [17] and a seam. Other 
numbers of PFs in MTs ranging from = 8 up to iV = 19 
have been observed as well [TTIH^ in the form of two- 
start or four-start helices corresponding to ft. = 1/N (a 
helical pitch of d) or ft, = 2/N (a helical pitch of 2d), 
respectively. In view of these different possible structures 
we want to study MT growth kinetics also as a function 
of the geometry parameter h [M] . 

In order to simulate the stochastic non-equilibrium 
growth dynamics of the model we use the Gillespie algo- 
rithm |23| . which implements a continuous time Markov 
process with the rates introduced above. For the simu- 
lations we use parameters, fcon — 200 [1/min], fcoff — 50 
[1/min], a monomer size d — 8 nm and a temperature 
ksT = 4.1pNnm corresponding to room temperature. We 
used a load distribution factor 6p = 1 and performed 
simulations both for energy distribution factors 9i = 
(interaction energy affects on-rate) and 9i = 1 (on-rate 
diffusion-limited). For each set of parameters, we average 
over 100 runs. 

Polymerization cycles, growth velocity, and stall 

force. — The polymerization kinetics is completely de- 
termined by the configuration of the N PF ends: the abso- 
lute length of each PF does not enter the rates ([2|) and we 
can neglect the possibility of vacancies or holes within a 
PF, which have not been observed experimentally. There- 
fore, the growth kinetics can be described by transitions 
between all possible states of the N PF ends. In general, 
each such state of the filament end can be described by a 



set n — (ni, ...,nN-i) oi N — 1 integer monomer number 
differences between neighboring PFs with Ui as monomer 
number difference between PFs i and i + 1. Transitions 
between these states happen by monomer addition and re- 
moval with rates ([2]). Monomer addition (removal) at PF 
i leads to changes Aui = +1 and An^-i = —1 (An^ = — 1 
and Ani_i — -1-1). The transition rates ka^^n^ between 
two states ni and n2 are only non-zero if both states are 
related by addition or removal of a single monomer and the 
corresponding rates are determined by ([2]), which depend 
on force and lateral interactions. 

During polymer growth, layers of monomers are added 
and eventually layers with N monomers are completed 
upon addition of a monomer. After addition of L complete 
layers the filament end attains the same configuration as 
initially. Therefore, each completion of L layers closes 
a polymerization cycle Cl of transitions in the network 
of states n. A non-zero average polymerization velocity 
V — (i) implies that layers are added with a non-zero 
rate and is equivalent to the existence of stationary cycle 
fluxes in the network of states n. We can calculate these 
stationary cycle fluxes using general theorems derived for 
the kinetics of chemical networks . 

We denote the stationary cycle flux for a polymeriza- 
tion cycle cl completing L layers by J{Cl). Likewise, 
the opposite cycle removing L layers is called and the 
corresponding stationary flux is J{C2)- For an arbitrary 
cycle cl — (ni, ^2, nM, f^M-i-i = ni) of length M and 
completing L layers {M > NL), the ratio of stationary 
cycle fluxes in forward and backward direction is given 
exactly by the ratio of products of transition rates along 
the edges of the cycles 

In addition, we can establish a general relation between 
the average growth velocity v and the stationary fluxes 
along cycles Cl'- The total stationary net flux along the 
fundamental set of all cycles Ci for completions of single 
layers gives the mean time to complete a single layer and, 
therefore, the mean growth velocity as 

V = dJ2{J{Ct) ^ J{C^)) (4) 

Cl 

Summation over single layer cycles Ci is sufflcient because 
they form a fundamental set |25| . i.e., linear combinations 
allow to represent cycles for an arbitrary number of L 
layers. 

At the stall force, the growth velocity v vanishes, i.e., 
J{Cl) — J{Ci), for all polymerization cycles for single 
layers. Because such cycles form a fundamental set, it fol- 
lows that all stationary net polymerization cycle currents 
have to vanish, i.e., J(C^) — J{C^) for all cycles Cl- 
This leads to the conclusion that the different fllament 
end states are in detailed balance at the stall force, as has 
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been conjectured in [TU]. According to the relation ^ we 
obtain the following Wegscheider condition 27 for any of 
the polymerization cycles Cl 



M 



IT. ^ 



= 1. 



(5) 



which will lead to an exact expression for the stall force. 

We first consider the cycle C^. Addition of exactly 
L layers requires M+ = (M -I- NL)/2 attachment and 
M_ = (M — NL)/2 detachment transitions, where M = 
M+ + M_ and NL = M+ - M_. Regardless of load 
and energy distribution factors, the thermodynamic con- 
straint ([T]) requires 



^fii,ni + i _ fcon ^±(ei\Aii,i+i\~F\Axi,i + i 
-^IjUi ^off 



D/ksT 



(6) 



for each attachment (+) and detachment (-) transition in 



C^. Therefore 



Jei\AiM\-F\AxM\)/kBT 



(7) 

where Mm = Y.t\ |A^m+iI " Efll = LNd is 

the total net gain in lateral contact length and IS.xm — 
Y^^i |Aa;i,i+i| — Y^^i |Aa;i,i+i| = Ld is the total net ad- 
vance of the leading PF. As a result, we find for the ratio of 
products of transition rates along the edges of the forward 
and backward cycles the simple result 



(fcon/fcoff) 



LN {eiLNd-FLd)/kBT 



(8) 



From the condition ([5]) we then obtain an exact expression 
for the stall force. 



install = N 



(9) 



Based on the assumption of detailed balance the same re- 
sult has also been obtained in ref. [13] . The stall force is a 
linear function of the number N of filaments in the ensem- 
ble and increases linearly with the lateral interaction. It is 
independent of all load or energy distribution factors of in- 
dividual polymerization steps. Moreover, the stall force is 
independent of the ensemble architecture. Our derivation 
shows that geometry independence not only means that 
the result ^ is independent of the parameter h for an ar- 
rangement with constant offset between neighboring PFs 
but that we obtain the same result © for the stall force 
for completely arbitrary arrangements of PFs relative to 
each other. Therefore, we also expect the stall force of a 
bundle of interacting actin filaments not to depend on the 
precise relative arrangement of actin filaments, which is 
hard to control in experiments. 

In the framework of polymerization cycles the so-called 
one-layer approximation introduced in 115| is equivalent 



to a "one-cycle" approximation, which restricts the sum 
in Q to a contribution from a single cycle dominating 
the sum in the limit of large forces close to the stall force. 
Therefore, the exact expression ^ for the stall force, at 
which all cycle currents become zero, is also recovered in 
the one- layer approximation in ref. |15) . This indicates 
how a systematic improvement of the one-layer approxi- 
mation could be achieved by inclusion of more polymer- 
ization cycles, which will leave the result ([SJ for the stall 
force unaffected. 

Experimentally, force-velocity curves have been mea- 
sured [3l[7] but the stall force is not directly accessible. 
Nevertheless, an exact result such as ^ can help to con- 
strain the analysis of experimental data. 

We conclude this section with a short discussion of the 
effects of GTP-hydrolysis on the staU force of MTs. GTP- 
tubulin attaches and can hydrolyze within the MT to 
GDP-tubulin, which gives rise to different off-rates fcoff,T 
and fcoff,D for GTP and GDP monomers, respectively. 
This leads to a coupling between polymerization cycles 
and hydrolysis, if the probability po that a monomer 
at the filament end is of GDP-type becomes non-zero. 
For small probabilities pd, effects from hydrolysis can 
be included approximately by using an effective off-rate 
fcoff,cfi = PTkoS,T +PDkoS,D in the result (O for the stall 
force, where pt = 1— is the probability that a monomer 
at the filament end is of GTP-type. Because hydrolysis 
destabilizes the filament and feoff, t < feoff, d, hydrolysis will 
generally reduce the stall force. For single actin filaments, 
the effect of hydrolysis on the force-velocity relation has 
been calculated in ref. [5S]. Hydrolysis and a non-zero pd 
also give rise to catastrophes as soon as the entire GTP- 
cap of a MT becomes hydrolyzed. 

Simulation results. — In simulations we can explore 
not only the stall force but the full force- velocity relation. 
We characterize the growth process by the average growth 
velocity v = (x). We determine the force- velocity relation 
as a function of both the PF interaction (per length) ei 
and the geometry parameter h. From the force-velocity 
curves, we determine the stall force numerically and inves- 
tigate how the stall force depends both on PF interactions 
and geometry parameter. Some of our simulation results 
have also been obtained in ref. [Tl] using a fixed time step 
Monte-Carlo algorithm. 

Force-velocity relation. We first describe results for 
the force- velocity relation of = 13 PFs corresponding to 
a MT using the rates ([2]). The shape of the force- velocity 
relation depends on the PF interaction ei [121 [TB]. For 

= all force-velocity curves end at the same velocity 
w(0) = (i(feon — feoff) for zero force, independently of the ge- 
ometry parameter h, see fig. [2^. Simulations confirm that 
force-velocity curves for different geometry parameters h 
exhibit the same staU force Fstaii = N{kBT/d) ln(fcon/feoff) 
as predicted in cq. The shape of the force- velocity 

curves between F — and the stall force F = install, how- 
ever, depends on the geometry parameter h, as shown in 
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Fig. 2: Geometry-dependence of the force-velocity relation for A'' = 13 PFs and two different values of the lateral interaction, 
(a) ei — OpN (curves are independent of the energy distribution factor ^i) and (b),(c) ei — S.OpN, ((b) 6i — 0, (c) 6i = 1). For 
each lateral interaction, different geometry parameters h (in units of nm and for monomer size d = 8nm) are compared. Insets 
are logarithmic plots. 



fig. [2^: for forces F > and small h, the kinetically limit- 
ing step is the insertion of the first monomer to an almost 
flat configuration. This rate-limiting step becomes faster 
for increasing h because the increase Aa; of the leading 
tip becomes smaller. This results in steeper force-velocity 
curves for decreasing values of the geometry parameter h, 
as can be seen in fig. [5^. 

For £; > 0, also the velocities at zero force v{0) depend 
on the geometry parameter h and on the energy distribu- 
tion factor 01. For increasing h, the zero-force velocity v{0) 
increases because the attractive PF interaction accelerates 
growth by reducing the off-rate (for 9i = 1, see fig. (He)) 
or increasing the on-rate (for 9i = 0, see fig. [IJb)). In- 
creasing the on-rate {9i = 0) leads to a much stronger 
effect. On the other hand, the simulations confirm that 
the stall force remains identical for different geometry pa- 
rameters h and for the different energy distribution factors 
6i in fig. ^h) and 0; = in fig. ^c) as predicted by eq. 
([9]). Nevertheless, it is not possible to simply conclude 
that force-velocity curves become increasingly steep for 
increasing geometry parameter h. Only for /i = 0, strong 
suppression of the rate- limiting first insertion step by force 
always results in the steepest force- velocity curves, see fig. 
[21 Then, the stall force is hard to determine because the 
force-velocity curve is extremely fiat with measured zero 
velocity over a range of higher forces. 

Stall force. It is possible to directly check the above 
exact analytical result (O for the stall force in simulations 
by the condition u (install) = 0, i.e., by determining the 
force where the average growth velocity in the simulation 
vanishes. In simulations this is done by applying a lin- 
ear interpolation to data points of force- velocity curves in 
vicinity of the stall force. 

In fig. [3] we show the simulation results for the stall 
force as a function of the lateral interaction ei for (a) N = 
13 and different values of the geometry parameter h and 
for (b) fixed h and different values of N. The simulation 
results clearly show a linear increase with ei and N and 
confirm the analytical result (|9]). In particular, we find no 
dependence on the geometry parameter h. 
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Fig. 3: Stall force (in pN) as a function of lateral interac- 
tion energy ei (in pN) for Oi — 1. Points represent simula- 
tion results, lines the analytical result ((9}. (a) For A'^ = 13 
PFs and different geometry parameters h = ^d/N ~ 0.934nm, 
h = d/N ~ 0.615nm, random h, h — l.Onm (from bottom to 
top). Curves are shifted by 0, 20, 40 and 60 nm respectively, 
(b) For a geometry parameter h — d/N and different PF num- 
bers A'^ = 4,5,6,8,10,11,13 and 15 (from bottom to top). 



The derivation of (jH]) was not limited to PF arrangement 
with a constant offset h between neighboring PFs but pre- 
dicts the same stall force for completely arbitrary arrange- 
ments of PFs. We checked this prediction by simulations 
of several random arrangements with random displace- 
ments between neighboring filaments, see fig. [3^ (green 
line). 

Dependence of velocity on the geometry parameter. 
The simulations allow us to explore how the growth veloc- 
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Fig. 4: Growth velocity as a function of the geometry param- 
eter h (in nm and for monomer size d — 8nm) for A'' = 4 (left) 
and TV = 13 (right) PFs for different values of the force F (in 
pN). The lateral interaction energy is ei = 3.0pN and the lat- 
eral energy distribution factor 6i = 1 for both top pictures and 
6i = for both bottom pictures. For vanishing force F — the 
curves exhibit A'^ — 1 maxima at ^ = i/N with i = 1, .., N — 1. 
With increasing force, some maxima vanish and new maxima 
can emerge. 



ity V depends on the geometry parameter h (0 < h < 1) 
for different load forces F, see fig. HI The curves v(h) are 
symmetric with respect to the axis h = 1/2. For e; > 
and zero load force F = 0, the curves v(h) exhibit — 1 
maxima corresponding to relative displacements h — i/N 
with i = 1, .., N — 1. At these values of h, the relative 
displacement h is commensurate with the monomer size d 
such that polymerization cycles are possible, where all sub- 
sequently attached monomers gain the same lateral inter- 
action energy. This avoids rate-limiting attachment steps 
and leads to optimal polymerization velocities. The max- 
ima are very pronounced for energy distribution factors 
01—0, where the on-rate is exponentially increased by 
lateral interactions and rather broad plateaus for 9i = I. 
This imphes that MT models using h=l/l3 [lOlfTT] over- 
estimate the growth velocity as compared to the actual 
three-start helix with h = 1.5/13. 

With increasing load F the height of the maxima de- 
creases, maxima can vanish or become minima, and new 
local maxima can appear. Upon approaching the stall 
force F w Fstaii, all curves v(h) become flat, which sup- 
ports the analytical result of an ^.-independent stall force. 

Conclusion. — Based on the concept of polymeriza- 
tion cycles we obtained the exact result (j9]) for the stall 
force of polymerizing ensembles of rigid protofilaments 
with lateral interactions. The stall force is a linear func- 
tion of the number N of filaments in the ensemble and in- 
creases linearly with the lateral interaction. On the other 
hand, the stall force is independent of the geometry of the 



ensemble and load or energy distribution factors. These 
results have been confirmed by simulations using the Gille- 
spie algorithm. Simulations also show that the shape of 
the force-velocity relation exhibits a pronounced depen- 
dence on the ensemble geometry below the stall force. Our 
results are relevant for the interpretation of experimental 
data on the force- velocity relation in microtubule polymer- 
ization and in the polymerization of bundles of interacting 
actin filaments or microtubules. 
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